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We study spectral properties of the Fokker-Planck operator that describes particles diffusing in 
a quenched random velocity field. This random operator is non-Hermitian and has eigenvalues 
occupying a finite area in the complex plane. We calculate the eigenvalue density and averaged 
one-particle Green's function, for weak disorder and dimension d> 2. We relate our results to the 
time-evolution of particle density, and compare them with numerical simulations. 
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In contrast to closed quantum systems, classical sys- 
tems often have dynamics generated by non-Hermitian 
operators. In this paper we develop general techniques to 
study the spectral properties of random non-Hermitian 
operators, and apply them to the Fokker-Planck (FP) 
operator that describes diffusion and advection of classi- 
cal particles in a spatially random but time-independent 
velocity field: 



dtn = C 



DV'^n - V ■ (Vn), 



(1) 



where n is the concentration of particles, D the molecu- 
lar diffusivity, and V the background velocity field. The 
non-Hermitian character, in this case, is due to the ad- 
vection term. 

The statistical behavior of the scalar field n obeying 
this FP equation has been investigated over a long history 
, in the context of turbulent diffusion and anomolous 
diffusion in random media. However, for general veloc- 
ity distributions, little appears to be known about spec- 
tral properties of the FP operator, such as the density 
of eigenvalues (DOS). An exception is the special case 
of potential flow, with V = V(j>. In this instance, there 
is a similarity transformation which maps the FP equa- 
tion exactly to a Schrodinger equation in imaginary time 
The eigenvalues of the FP operator with potential 
flow are therefore real and negative. In one dimension, 
it is possible to express any velocity field in terms of 
a potential and to transform the FP equation in this 
way. Moreover, anomolous diffusion in one-dimensional 
systems with random flow has been shown to be con- 
nected with logarithmic singularities of the DOS and of 
the eigenstate localization length as the eigenvalue ap- 
proaches zero In two dimensions, in the opposite 
case of incompressible flow (V.V — 0), for which there is 
no similarity transformation to a Hermitian operator, lo- 
calization of eigenfunctions of the FP operator has been 
studied numerically Q . In addition, an analogy between 
the classical FP equation and the quantum random flux 
problem has been analyzed 

In general, the eigenvalues of a non-Hermitian FP op- 



erator occupy a finite area in the complex plane, rather 
than being restricted to the real-axis. This fact, despite 
the similarities in other respects between the FP and the 
Schrodinger equations, renders inapplicable ||^ the stan- 
dard perturbation expansion of Green's functions, used 
for disordered quantum problems. Furthermore, saddle- 
point techniques developed for non-Hermitian ran- 
dom matrix ensembles are too specialized to be 
appropriate for random FP operators. 

Very recently there has been considerable interest in 
properties of random non-Hermitian operators [p^|-p^, 
with a range of motivations, including the study of open 
quantum systems and the motion of flux lines in su- 
perconductors [|l2|. Against this background, a better 
understanding of their spectral properties and of calcu- 
lational methods is clearly desirable. 

In this paper we describe a general scheme, based 
on a diagrammatic expansion, to compute the disorder- 
averaged Green's functions and the DOS of random non- 
Hermitian operators. Similar ideas have been proposed 
in the context of random matrix theory by Feinberg and 
Zee p^. We apply the technique to the FP operator of 
Eq. (|l^ calculating the shape of the support of the DOS, 
and the eigenvalue density itself. We also compare our 
analytical results with numerical calculations. 

The particular FP operator we consider has constant 
diffusivity D and a quenched random velocity field V(x). 
Note that this is the opposite limit to that in the model 
recently discussed by Kraichnan Q|, which has infinitely 
short time correlation in velocities. Time-independent 
fiows can be established in physical systems such as 
porous media [||. 

We take the velocity field to be gaussian distributed, 
with zero mean, and variance 



(v„(k)v^(k')) = ri(<5„^ - -^)S{k + k' 



r2(^)5(k + k') 



(2) 



where V(k) = (27r) J d'^xe **'-^V(x), angular brackets 
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denote the ensemble average, and Fi and r2 represent the 
strengths of the transverse and longitudinal parts of the 
velocity field. We take the spectrum of velocity fluctua- 
tions to have a short- wavelength cut-off, A, and consider 
a system of volume il. The special case of mixed flow, 
Fi = r2 = F, leads to substantial simplifications: for 
clarity of presentation, we describe calculations only in 
this limit, but state results for the general problem. 

At a complex frequency cj, the dimensionless combi- 
nation 7(w) = (F/|wn(|tj|/D)''/2+i is a measure of the 
disorder strength. The fact that, as w ^ 0, 7(0;) ^ for 
d > 2 and 7(0;) — * 00 for d < 2, identifies d = 2 as the 
upper critical dimension [||. 

Our aim is to study spectral properties of the FP op- 
erator, Cfp. We do so via the ensemble-averaged Green's 
function 



1 



c 



fp 



(3) 



Let (La I and \R\) be left and right eigenvectors of Cfp 
with eigenvalue A, and let |p) denote a plane- wave basis 
state with wavevector p. The ensemble-averaged spectral 
density 



C{p,Lu) = {Y^{p\R,){L,\p)S{u;-\)) 



(4) 



is diagonal in this basis, because averaging restores trans- 
lational invariance. From it we can obtain the time 
Green's function, or particle density, (n(r,i)), evolving 
according to Eq (|l]) with initial condition n(r, t — 0) — 
Sir): 



(n(r,t)) = (2^)- 



ip.r 



*C{p,Lo). (5) 



In the same basis, the diagonal elements of g{uj), which 
we compute, are related to C(p,u;) by 



LU — X 



(6) 



Analytic properties of the Green's function depend on 
the eigenvalue density, p(A), in the complex A-plane. In 
particular. 



(7) 



Thus, in the complex cj-plane, g{Lu) is non-analytic ev- 
erywhere that the eigenvalue density is non-zero. Stan- 
dard techniques for calculating disorder-averaged Green's 
functions via perturbation theory yield only the part of 
g{uj) that is analytic in uj, together with its analytic con- 
tinuation inside the support of Q. In the special 
case of pure potential flow, (Fi — 0), the eigenvalues of 
Cfp all lie on the negative real axis and g{uj) is analytic 
elsewhere. Under these circumstances one can compute. 



for example, p{oj), in the usual way, from the discontinu- 
ity in g{Lu) across the real axis. By contrast, for mixed 
flow (as we shall show), the eigenvalues of Cfp fill a fi- 
nite area in the complex plane, the analytic part of g{uj) 
contains limited information, and a new approach is re- 
quired. 

To this end, for a general C, we construct an matrix Ti. 
which: (i) has twice the dimension of C; (ii) is Hermitian; 
and (iii) has an inverse that contains g{u}) as one of four 
blocks. Specifically, with A = uj — C, 



A 



A^ 



(8) 



The inverse, G = TL ^, exists for real e > and is 



G 



Gil G12 
G21 G22 

{e + AA^y 
A\e + AA^) 



A{e + A^A)-^ 
-e(e + At^)-i 



(9) 



Since Ti is Hermitian, G can be calculated using estab- 
lished methods, and from it we can obtain 



g{uj) = lim(G2i). 



(10) 



This approach parallels recent work by Feinberg and 
Zee ||l^ and is somewhat different from that taken 
in other calculations on spectral properties of non- 
Hermitian operators: in the present notation, Sommers 
and collaborators ||,|l^ focus on det(G'ii), while Efetov 
separates C into Hermitian and anti-Hermitian parts. 

A diagrammatic expansion for (G), and hence g(ti^), 
follows from writing 



where 
Ho 







Ti. — Ti-o + Ti-i , 



and Til 



A 
A^ 



(11) 



(12) 



The series for G, in powers of Hi and of G^"^ = Hq ^, 
involves two propagators, Gf^ and G22^ , and two vertices, 
A and A^ . As usual, it is convenient to introduce a self- 
energy, S, and proceed via the Dyson equation. 



(G) = Go + GoS(G>, 



where 



Go 



1/e 

- 1 



(13) 



(14) 



To illustrate the approach, we consider first the asym- 
metric gaussian random matrix ensemble of Ref Q : the 
N X N real matrix [J] has the distribution P[J] cx 
exp{-NTr{Jj'^ - t J J)/[2(l - r^)]) so that non-zero co- 
variances are 
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< Jfc. >=l/N, < J,kJk, >= t/N. 



(15) 



The fully asymmetric problem, in which t — 1 and Jik 
and Jfci are statistically independent, was first studied 
by Ginibre and has been treated using a Green's 
function methods in Ref [ p^ . 

Setting C — J, the leading contribution to the self 
energy at large N comes from the diagrams of the self- 
consistent Born approximation (SCBA), shown in Fig|^. 







FIG. 1. Self-energy diagrams in the SCBA. Single and 
double lines denote the two components of the bare Go: 1/e 
and —1, respectively, while the internal digrams describe the 
four full Green's functions. 



From these wc obtain 

S = 



(G22) r{G2i)-Lu 
t{Gi2)-^* (Gn) 



(16) 



Solving eq.(jT^), we recover the results of Ref In 
particular, inside the support of the p(ti^), defined by 
[x/{l + r)]2 + [y/{l - r)]2 < 1, we find 



N-^Trg{uo) = ( 



ly 



l + T 1 - r 



), uj = x + iy, (17) 



and hence a constant DOS, p(cj) = [it{1 — t^)] ^ . 

With this preparation, we return to the Fokker-Planck 
operator of Eqs (0) and (|), taking Fi = Fj = F. We find, 
at weak disorder and for dimension d > 2, that the self- 
energy is again given by the SCBA. Corrections are small 
by a factor of the dimensionless disorder strength, ^{x); 
neglecting these, we obtain the self-energy (diagonal in 
wavevector, p) 



(18) 



where cUp ^ lu + Dp^ . 

Before discussing the self-consistent solution to Eqs 
( |l3| ) and (18), we note that one can arrive at the same 
point by expressing G in terms of a functional integral, 
averaging using replicas, and making a simple decoupling 
approximation. In view of the 2x2 block structure of G, 
it is natural to introduce two complex fields, 771 and 772, 
with rj^ = {ri\, rjl), and write for j,k = 1 or 2 



G 



jk 



i(Z- 



(19) 



where the normalisation and action are 

Z = J V[ri]e-'^ and S = r/t^?] 



(20) 



The average over the velocity field, denoted by (. . .) gen- 
erates terms in the action of the form ri\ Q,ryi,/3??2 /3'?2,a, 
where a and (3 label replicas. Approximating these by 
setting 



we arrive at 



(21) 



(22) 



with 5* [G"] ^ — S, where E is again given by Eq |T^. 

Proceeding now to the evaluation of the self-energy, we 
find 



(Gfi 



(Gfi) 



Up I ~r -LJp 



\ujp\^ + Bp 



pi ' 

' . (G^2) , (23) 



\ujp\^+Bp 



where Bp = {e + Fi){\ + F2) and 



Fi = -F r d\{Gl^), F2 =F t d\q'{G\,). 
Jo Jo 

Analysing this system of equations in the limit e — > 0, we 
obtain 



\ujp\^ + bp^ 



(24) 



where the value of b is determined via behaviour of the 
integral 



Jo 



\ujq\^ + bq^ 



(25) 



If I{b = 0) < 1, then b — 0; otherwise b is the (real and 
positive) solution to the equation I{b) — 1. The former is 
the case if cu is not close to the negative real axis. In that 
event, with & = 0, 5^(w) = uj~-^ and dg^{uj)/dui* = 0, so 
that p{llj) — 0. Alternatively, if w is close to the negative 
real axis, 6 > 0, ^^(w) is not analytic in to and p{uj) ^ 0. 

Thus the boundary to the support of the density of 
states, 2/b(x), satisfies the equation 



.A 

7(6 = 0) EE F / d'^q- 
Jo i 



\x + Dq^y+yl 
For small F, d > 2 and x < 0, we find 

yB(x) = ±'^^{x)\x\^±V\x\'"^ 



= 1- (26) 



(27) 
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where Sd is the surface area of a d-dimensional unit 
sphere. The DOS in the region a; < 0, \y\ < ys is 

D 



(27r)d+ir|a;| 



(28) 



and elsewhere is zero. Thus the eigenvalues occupy 
a wedge-shaped region in the complex plane, centered 
around the negative real axis. The x-dependent width 
of this region can be understood simply: if one assumes 
that it is proportional to T, dimensional analysis implies 
ysix) oc 7(x)|x|. Similarly, the ^-dependence of p{uj) 
follows from the requirement that J dyp{x -t- iy) take the 
value it has in the disorder-free system. 

Turning to the time-evolution of the particle density, 
we find that the time Green's funtion, defined in Eq (|^), 
is the Fourier transform of a product of two factors: 



{n{r,t)) 



1 



(29) 



where ujq — yB{x), evaluated at x = —Dp^: cjo = 
{'KSdTq'^)/{2D). The first factor, e~°P^\ is the famil- 
iar consequence of simple diffusion; the second factor, 
s,m.{ujQt)/<jjQt, arises from advection. At weak disorder, 
when the SCBA is a controlled approximation, the ad- 
vective factor differs significantly from 1 only where the 
diffusive factor is small, so that (r^) ~ Dt. By contrast, 
at strong disorder (when the SCBA is simply a mean- 
field approximation), it is the advective factor that sets 
the width of of the density profile at short times, and 
(r^) {Tt/ D)'^/'^ . This time-dependence arises because 
advection and diffusion contribute equally to particle mo- 
tion in this regime. 

For general random flow, with Fi ^ F2, an extension 
of our approach again yields Eqs (|2^), ( p8| ) and (29), but 
with Fi replacing F. Thus, in particular, for potential 
flow (Fi — 0) we reproduce correctly the fact that all 
eigenvalues are real. 

To test the theory developed above, we have calculated 
numerically the eigenvalues of the FP operator for mixed 
flow, discretized on a square lattice. Theory (adapted 
to the discretized FP operator) and simulation are com- 
pared in Fig (H), for 50 realisations of a 32 x 32 lattice 
with D = 1 and F = 0.25. 

The calculated shape of the support of p(w) matches the 
data well; the fact that a finite fraction of the eigenval- 
ues have vanishing imaginary part is a finite-size effect, 
analyzed in detail in Ref jlj]. 
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FIG. 2. Distribution of eigenvalues (dots) in the complex 
plane, and calculated boundary (o) of p{io). 
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